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Abstract 

We discuss physical implications of the explicit method in numerical anal- 
ysis. Numerical methods have there own condition for causality, known as the 
Courant-Friedrichs-Lewy condition. It is proposed that numerical causality 
merges with physical causality as the grid interval size approaches zero. We 
discuss the implications of this proposition on the numerical analysis of the 
wave equation. We also show that, insisting on physical causality, the numeri- 
cal analysis of Schrodinger's equation implies that the minimum space interval 
should satisfy Ax > oqAc, where Ac is the reduced Compton wavelength and 
ao is a constant of the order unity. 
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1 Introduction 



Partial differential equations (PDEs) are ubiquitous in all fields of science and are 
used to model physical systems. In numerous cases, obtaining the exact solutions 
of these equations is not possible and one has to resort to numerical approximations 
in order to solve these systems. Numerical methods therefore have a wide range of 
applications in numerous scientific fields. In physics, for instance, these methods are 
used to obtain approximate solutions of important differential equations such as the 
diffusion and wave equation. 

In this article we attempt to demonstrate how a formalism used to solve PDEs 
in applied mathematics can be employed to further our understanding of known 
fundamental ideas in physics. The numerical solutions of widely used equations 
in physics, such as the wave and heat equation, are known and there are different 
methods to solve them. We shall use the explicit finite difference method to obtain 
solutions of the wave and Schrodinger equations. We will use the Von Neumann 
stability analysis to discuss the limits on the grid intervals. This analysis typically 
leads to limits on the space and time intervals of the numerical grid in order to ensure 
stability. We will discuss these limits while emphasizing on physical causality. 

The notion of causality is one of the foundational ideas in modern physics. In the 
context of Einstein's special theory of relativity it implies a Universal limit on the 
speed with which cause and effect are related. This universal limit being the speed 
of light. It is the limiting speed with which information can travel in the Universe. 
In other words, causality implies that the effect of an event in space-time should 
belong to its future light cone. The stability requirement of the explicit method in 
numerical analysis typically leads to a condition of causality, known as the Courant- 
Friedrichs-Lewy (CFL) condition. We propose in this article that in the continuum 
limit numerical causality unifies with physical causality. We will show that insisting 
on physical causality in the explicit method of the Schrodinger equation can lead to 
limits on the minimum grid interval sizes which are essentially consistent with those 
obtained in quantum mechanics. 

The paper is organized as follows: In section 2 we briefly introduce the finite 
difference methods used in numerical analysis to solve partial differential equations. 
Section 3 includes a discussion on the wave equation and the stability condition 
which implies an upper limit on the speed of waves. In section 4 we show that the 
convergence and causality in the explicit method of the Schrodinger equation imply 
the minimum grid interval size to be of the order of Compton wavelength. Section 5 
includes a general discussion on the properties of the explicit method that accord it 
a causal structure. We conclude in section 6. 
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2 Numerical analysis and discretization 

Finite difference methods (FDMs) are frequently used in numerical analysis to solve 
partial differential equations. In these methods a PDE is solved on a grid of discrete 
points with, for instance, Xn = riAx and tj = jAt, where n and j are integers [1, 2, 3]. 
Continuous derivatives are then replaced by finite difference approximations. The 
solution of complex differential equations can be obtained by dividing space-time into 
a discrete set of points in this manner. The finite difference equation, for example, 
can be a forward, backward or centered difference. Each choice needs to be tested in 
various ways such as its stability, accuracy and speed. 

A numerical method is stable if the errors in the method do not increase as 
the solution moves forward in time. There can be several sources of errors in a 
numerical method. These include, for example, the truncation and round-off errors. 
The round-off error results from finite precision representation of real numbers on 
a computer. Truncation errors arise from representing continuous derivatives with 
finite differences which involves truncating terms from the Taylor series. The leading 
terms in the truncation error determines the accuracy of a FDM. 

What makes a particular numerical method feasible is its consistency, i.e., whether 
the finite difference equation (FDE) approaches the PDE as the interval size vanishes. 
In other words, a numerical scheme is consistent if the truncation error goes to zero 
as Ax, At — )■ 0. Similarly, convergence is also important and requires the solution of 
the FDE to converge to the exact solution in the limit Ax, At — >■ 0. According to the 
Lax equivalence theorem, if a finite difference method is consistent (FDE approaches 
the PDE as grid interval size approaches zero) and stable (errors do not grow) than 
convergence (discrete solution converges to the actual solution as grid interval size 
approaches zero) is guaranteed. In short, consistency and stability imply convergence. 

The most convenient way to test the stability of a finite difference method is 
by employing the Von Neumann stability analysis. In this method we test how a 
Fourier mode behaves on the grid. Compared to other methods, such as the matrix 
method, this analysis is a convenient way to test the stability of a method. One short 
coming of this method is that it ignores the boundary conditions. This is why the 
limits obtained are necessary but not sufficient to guarantee stability. But this short 
coming, in a way, is useful in our case since the limits we get would be independent 
of the boundary conditions of the problem. In FDMs, a first order derivative in space 
and time can be approximated by a forward time and space difference as follows 

du Uj+i — Uj 
dx Ax 
du u^'^^ — 

'di ~ At 

where i = 0, 1,2 ... M and n = 0, 1,2 ... N. Amongst several methods used, two 
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Figure 1: Stencil for the FDM of the 1-D wave equation. It represents the expression 
given in equation (9). 

important ones are the exphcit and imphcit methods. Each of these has its own 
benefits in terms of stabihty, speed and accuracy of the solution. The explicit method, 
which is a forward difference method, is fast but less stable since it usually requires 
limits on the space and time step sizes. In the explicit method the state of the system 
at a time is evaluated from the state at time t„ or a prior time (tn-i, ■■), namely 



The implicit method, which is a backward difference method, is slow and more 
numerically stable. In this method the state of the system at a later time is not an 
explicit function of the state at earlier times. For this case, therefore, the state of the 
system at a later time would typically be given by 



The factor that makes this method more stable is that the time step is usually not 
constrained by a stability condition. A larger time step At can therefore be used 
to solve the system. However the above equation is usually a non-linear system of 
equations and needs additional time to solve. Therefore this method is slower than 
the explicit method. We will focus on the explicit method in the following sections. 

3 Wave equation: Speed Limit 

In this section, we show that employing the explicit finite difference method in solv- 
ing the wave equation can lead to a limit on the propagation speed of waves. The 



f i^ni tn) ■ 



(3) 



Un+l = f{u 715 Un+l, tn+l)- 



(4) 



4 



wave equation is a hyperbolic PDE and describes the propagation of waves. Hyper- 
bohc PDEs describe systems with a finite speed of propagation. Consider the one 
dimensional wave equation: 

We can solve the wave equation using the explicit method by replacing the derivatives 
by centered differences in space and time (leapfrog method) as follows 

A^^ + ^^^^^)' 

where = u{xj,tn)- The truncation error in this case is O(Ax^) and O(At^) and 
arises from truncating terms from the Taylor expansion. The leading terms in the 
truncation error of this method are given by 

T," = ^v\AxY u,,,, - 4(At)2 utttt. (8) 
Dropping the truncation error we get the following approximation 



u]+^ = r- 



K+i - « --i) + 2(1 - r>," - \ (9) 



where r = vAt/ Ax. To study the stability of this method we use the Von Neumann 
or Fourier analysis. We look for solutions of the form 



^ ^^K3^-)G^{At, Ax, k), (10) 



and get the following equation for the growth factor 

- 2[1 ~ 2r^ sm\kAx)]G + 1 = 0. (11) 
The requirement that IG] < 1 leads to the following CFL condition for stability 

Ax , , 

In 3 dimensions, the CFL condition for stability is given by [2] 



1 

- At 
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Aa;2 Ay^ Az^ 



-1/2 

= AVgrid, (13) 
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which for a cubic grid {Ax = Ay = Az) will be 



V < 



1 Ax 
7^ At 



AVgrid- 



(14) 



Here Avgrid is the average grid velocity and can be viewed as the speed with which 
information moves across the grid [4]. We can view equation (12) as a limit on the 
velocity of the wave. It requires the domain of dependence of the wave equation to 
be within the numerical domain of dependence. The question of making Ax and At 
smaller is an important one in numerical analysis. As At is made smaller the time 
for the computation increases rapidly. So, practicality also dictates At not to be too 
small. 

As described in section 2, in order for the explicit method to be feasible it should 
also be convergent. The truncation error for this method given in equation (8) is 
0(Ax^,At^) and vanishes in the limit Ax, At — )■ 0. However, as the grid interval 
size approaches zero we need to specify the relationship between Ax and At. We 
approach Ax, At — > such that the ratio Ax/At is kept constant. This is referred 
to as the refinement path [3]. 

Furthermore, the speed with which numerical information for the explicit method 
flows on the grid is different from the physical speed of information transfer. Herein we 
will assume that the numerical speed of information becomes equal to the physical 
speed of information as the grid approaches the continuum limit. In other words, 
numerical causality merges with physical causality in the continuum limit. Therefore, 
taking the limit Ax, At — > while keeping the instantaneous velocity Vgrid constant 
we get 



where the constant c is the speed of light. It is the limiting speed with which infor- 
mation travels in the continuum limit. Therefore, by assuming the numerical speed 
to be the same as the physical speed limit as Ax, At — )■ 0, we get a limit on the speed 
of waves which obeys causality. This limit also ensures the stability of the explicit 
method. 

Note that the numerical solution of Maxwell's equation is obtained using the Finite 
Domain Time Difference (FDTD) method and results in a conditional stability limit 
similar to equation (13) (see, for example, chapter 4 of ref. [5]). For the case of 
Maxwell's equation the condition (13) becomes an equality in the continuum limit. 

Causality: The domain of dependence (DoD) and range of influence (Rol) of a 
hyperbolic PDE is shown in Fig. 2. The DoD of a point in the solution domain is the 
set of points on which the solution at a particular point depends. Similarly, the Rol 
is the solution domain which is affected by the solution at a particular point. The 



V < Vgrid = C, 
=^ V < C, 



(15) 
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Figure 2: Figure shows the domain of dependence and range of influence for the 1-D 
wave equation. The gray region bounded by the sohd hne represents the continuous 
PDF whereas the dashed hne shows the boundary of the numerical solution. The 
CFL condition in equation (12) requires the DoD of the PDF to be contained within 
the DoD of the numerical solution. As the grid interval size approaches zero the 
slope of the dashed boundary representing the numerical PDF should approach the 
limiting speed with which information can be transferred across the grid. So in the 
continuum limit, f < c, if we insist that numerical causality merges with physical 
causality. 

CFL condition in equation (12) also ensures that the PDF domain of dependence 
stays within the numerical domain of dependence as shown if Fig. 2. The gray region 
is the DoD of the PDF which is contained in the DoD of the numerical solution. 
The numerical value of the solution at a node (xj,t„) depends on the values at the 
nodes lying within the numerical domain of dependence. The wave equation therefore 
obeys causality and, particularly, in the continuum limit, equation (15) guarantees 
that physical causality is always satisfied. 
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4 Schrodinger equation: Minimum Space and time 
intervals 



In this section, we show that insisting on physical causahty and convergence in the 
exphcit finite difference method to solve the Schrodinger equation yields a limit on 
the minimum grid interval size. The Schrodinger equation can be obtained from 
the non-relativistic limit of the relativistic expression for energy, E"^ = •p'^ + m?, by 
replacing relevant variables by operators as, — m ~ if = ihd/dt and p = h/i d/dx. 
The time dependent Schrodinger equation for a free particle in one dimension is given 

by 

d%jj ih d'^ip 

'dt^2^d^' ^^^^ 

The Schrodinger equation is similar to the diffusion equation but with an imagi- 
nary diffusion coefficient. This changes the properties of the solutions completely. 
Schrodinger's equation allows plane wave solution whereas the solutions of the diffu- 
sion equation are diffusive. We first employ the finite difference form of the Schrodinger 
equation (16) with forward difference in time and centered difference in space (FTCS) 
and show that it is unstable. Therefore, equation (16) can be approximated as 

*«'-^*? + ^^ + 0(A,A.^). (17) 



which yields 



At 2m (Ax) 



where, (3 = h/m and a = At/2{Ax)'^. The leading terms in the truncation error of 
this method are given by 

= \At - ^^i^^y (19) 

Consistency requires this error to vanish as the grid approaches the continuum limit. 
To find the growth factor G = G{At, Ax, k) we use the Von Neumann analysis by 
substituting the following expression in equation (18) 

^ e^'^(^'^^)G"(At, Ax, k), (20) 

and attain the following equation for the one step growth factor G 

G = l + 2ia(3[cos{kAx) - 1]. (21) 
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The condition for stability is that IG] < 1 should be true for all kAx. The magnitude 
of the growth factor in equation (21) is 

= 1 + [2(3a{cos{kAx) - 1)]^ > 1. (22) 

Therefore, the growth factor is greater than 1 for all kAx ^ and this method 
is unstable. The FTCS method along with some other explicit schemes to solve 
the Schrodinger's equation are known to be unstable [6]. For this reason the im- 
plicit methods are typically employed to solve it. The fully implicit method for the 
Schrodinger equation involves the centered difference in space, i.e. the right hand 
side of equation (17), at time tn+i- In the Crank Nicholson method the average of 
the explicit and implicit schemes is taken and this method is known to be uncon- 
ditionally stable for the Schrodinger's equation [7]. Although stable, this method is 
computationally expensive as it involves inversion of large matrices. Moreover, it has 
infinite speed of propagation of numerical information and is therefore not suitable 
if we insist on physical causality in the continuum limit. As mentioned earlier, this 
is the reason we will focus on the explicit method. 

One of the reasons for the instability of the FTCS method is that it is not cen- 
tered in time [8]. To overcome this an explicit finite difference scheme was proposed 
in [9] which is centered in time, unitary and conditionally stable. The following 
approximation is taken in this method 

2At 2m (Aa;)2 ' ^ ^ 

The stencil for this equation is similar to the one for the wave equation shown in 
Fig. 1. The truncation error for this case is 0(At^, Ax^). The leading terms in the 
truncation error of this method are given by 

T," = ^{Atf iJut - ^^(^^)' (24) 

The Von Neumann stability analysis leads to the following equation for the growth 
factor 

+ Aia(3[l - cos{kAx)]G -1 = 0. (25) 
The growth factors are therefore 



^1,2 = -2ial3[l - cos{kAx)] ±^/l- 4a2/32[i - cos{kAx)]'^. (26) 
The requirement |G| < 1 leads to the following stability condition [9, 10] 

a/3 < ^, (27) 
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which imphes 



At ~ m 



{2i 



The stabihty condition in 3 dimensions leads to the following limit [9] 
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For a cubic grid Ax = Ay = Az we get 

Ax2 2h 

> — . 

3At - m 

As before, defining Avgrid = Ax/ a/3 At, equation (30) becomes 

2h 



(29) 



{AVgridYAt > 



m 



or 



AEg„dAt > h. 



(30) 



(31) 



(32) 



The above equation implies that in order for the numerical solution of the Schrodinger 
equation be stable the grid interval must satisfy the above limit. The above limit 
is valid for any grid size but for the method to be viable this should particularly be 
true for the limit Ax, At — 0. As mentioned earlier, convergence is an important 
condition for a FDM to be viable. Therefore, in the limit Ax, At — )■ we approach the 
quantum realm and in addition make sure that our difference method is convergent. 

Refinement path: As we approach the limit Ax, At — we should also specify 
the relationship between Ax and At. As we approach the limit Ax, At — )■ 0, the 
truncation error approaches zero, and the grid velocity Vgrid remains constant. Here 
Vgrid is the speed with which numerical information flows on the grid. Assuming 
numerical causality approaches physical causality in the continuum limit this constant 
should be the speed of light. As the mesh interval size approaches zero this is the 
limiting speed with which information can be transmitted on the grid. Therefore, Ax 
and At — )■ such that Vgrid = c and the limits on the length and time intervals are 



At > 



2% 



(33) 



or 



Ax > 2v^Ac 



(34) 
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where Ac = h/mc is the reduced Compton wavelength. The above condition imphes 
that causahty requires a hmit on the minimum spatial grid interval size to be of the 
order or greater than the Compton wavelength. Another implication of these limits 
is that the truncation error in the FDE is 0((Ax)^) > 0{\\). Note that we have 
assumed a cubic grid to derive this limit and similar limits would be true for Ay and 

In reference [11], the stability condition for the Schrodinger equation is derived 
using the FDTD method. The limit in equation (28) for that case is h/m instead 
of 2H/m. Similarly, other higher order methods [12] also lead to stability condition 
similar to equation (28). Therefore the stability requirements on the limit on grid 
interval size for various schemes can be written as 

Ax > a^Xc, (35) 

where Oq is a constant of order unity and varies from one explicit scheme to another. 
Hence, in addition to viewing Compton wavelength as a limitation on measuring the 
position of a particle it can also be viewed as the minimum grid length that can 
be resolved in order to ensure causality in the stability analysis of the Schrodinger 
equation. Furthermore, assuming a minimum length in this manner renders the 
theory the same causal structure as the wave equation. 

5 The explicit method and its imphcations 

We have seen in previous sections that the explicit method can shed light on some 
important ideas in physics. It is important to note the difference between the explicit 
and implicit method. The explicit method is fast but less stable and requires small 
time steps for more accuracy. The implicit method is slower, more stable and allows a 
larger time-step. In addition to these factors, the explicit method has some important 
features that makes it different from implicit methods. Following can be the two 
reasons as to why the explicit method yields these important relationships in contrast 
to the implicit methods: 

1. The solutions at present time are explicit functions of solutions of past. This 
implies temporal causality. 

2. The solution at a grid location do not instantaneously effect the solutions at 
the entire grid [4]. This implies a finite speed of information transfer which is 
called the grid speed. 

We can better understand these points by considering the example of the heat 
equation which does not limit the speed of heat transmission. However, when using 
the explicit finite difference method in seeking the numerical solution of a differential 
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equation this might not necessarily hold. The heat equation, for instance, has the 
instant messaging property but it does not hold when we employ the explicit method 
to solve it [4]. The speed with which information travels through the grid is Ax/ At 
and not infinite. This, of course, implies that the explicit method is not the perfect 
method to solve the heat equation but at the same time it does render the system 
a causal structure. The implicit method, however, does have the instant messaging 
property and this is apparently the reason it is unconditionally stable. Therefore, 
the numerical propagation speed of information for explicit methods is finite whereas 
that of implicit method is infinite. 

The Schrodinger's equation is similar to the heat equation but with an imaginary 
diffusion coefficient. It is apparently a parabolic PDE but the solution it permits has 
properties of the solutions of hyperbolic PDEs. In other words although it satisfies 
the mathematical test of a parabolic equation which describes diffusive processes 
it allows plane wave solutions. Parabolic PDEs allow infinite speed of information 
transfer whereas this speed is finite for hyperbolic PDEs. In both cases, however, 
the explicit finite difference method when employed allows only a finite speed of 
propagation. 

6 Conclusion 

We discussed how the techniques used in numerical analysis might be more than tools 
to obtain solutions of partial differential equations. Finite difference methods and, 
in particular, explicit methods might help us further understand some fundamental 
ideas in physics. Explicit methods used in solving finite difference equations are often 
conditionally stable. The explicit finite difference method for the wave equation can 
lead to a limit on the speed of waves to be less than the grid speed. As the grid 
approaches the continuum limit, we assumed that the idea of numerical causality 
merges with physical causality and the speed of numerical information transfer be- 
comes equal to the physical speed limit, i.e., the speed of light. This leads to a limit 
on the speed of waves that ensures stability of the explicit method and also ensures 
physical causality. 

Similarly, the explicit method of the Schrodinger equation yields a limit on the 
minimum length and time interval. To ensure convergence of the explicit method, 
we approach zero interval size along the refinement path Vgrid = c. If we insist on 
physical causality in the continuum limit this speed should be the limiting speed of 
information travel, i.e., the speed of light. We thereby showed that the limits on the 
minimum spatial and temporal grid lengths are Ax > 2\/?)h/mc and At > h/mc^. 
This can be understood as follows. In the continuum limit, the grid interval, being 
a measure of the position of the particle, suffers from similar limitations as those 
obtained in quantum mechanics. These limits result from the Von Neumann stability 



12 



analysis. Therefore, the exphcit method endows a causal structure to the Schrodinger 
equation which is similar to that of the wave equation and as a result implies the 
minimum grid interval size to be of the order of Compton's wavelength to guarantee 
stability. 

Finally, we also discussed the importance of the explicit finite difference method in 
numerical analysis. The Von Neumann stability analysis used to test the stability of 
explicit methods implies that our conclusions are essentially independent of boundary 
conditions. The fact that this method preserves causahty and has a finite speed of 
information transfer can be two important reasons it yields important limits on the 
spatial and temporal grid interval sizes. Therefore, in our analysis we showed that 
insisting on physical causality and convergence of a finite difference method can lead 
to important ideas in fundamental physics. 
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